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Abstract: In the Next-to-Minimal Supersymmetric Standard Model (NMSSM), all 

singlet-dominated particles including one neutralino, one CP-odd Higgs boson and one 
CP-even Higgs boson can be simultaneously lighter than about 100 GeV. Consequently, 
dark matter (DM) in the NMSSM can annihilate into multiple final states to explain the 
galactic center gamma-ray excess (GCE). In this work we take into account the foreground 
and background uncertainties for the GCE and investigate these explanations. We carry 
out a sophisticated scan over the NMSSM parameter space by considering various experi¬ 
mental constraints such as the Higgs data, H-physics observables, DM relic density, LUX 
experiment and the dSphs constraints. Then for each surviving parameter point we per¬ 
form a fit to the GCE spectrum by using the correlation matrix that incorporates both 
the statistical and systematic uncertainties of the measured excess. After examining the 
properties of the obtained GCE solutions, we conclude that the GCE can be well explained 
by the pure annihilations bb and XiXi Ai Hi with A\ being the lighter singlet- 

dominated CP-odd Higgs boson and Hj denoting the singlet-dominated CP-even Higgs 
boson or SM-like Higgs boson, and it can also be explained by the mixed annihilation 
XiXi W + W Among these annihilation channels, XuXi ~> A\Hi can provide 

the best interpretation with the corresponding p- value reaching 0.55. We also discuss to 
what extent the future DM direct detection experiments can explore the GCE solutions 
and conclude that the XENON-1T experiment is very promising in testing nearly all the 
solutions. 
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1 Introduction 


The compelling evidences for the existence of Dark Matter (DM) from various cosmological 
and astrophysical observations have provided us a good portal in the search for new physics 
beyond the Standard Model (SM). One possible method to explore DM in the present Uni¬ 
verse is the indirect detection, which looks for the particles produced when DM annihilates 
in the DM halo. These particles include photons, antiparticles and neutrinos, and among 
them gamma rays have often been defined as the golden channel for DM indirect detec¬ 
tion since the signal can be traced back to the source. The Large Area Telescope (LAT) 
onboard the Fermi Gamma-ray Space Telescope, due to its unprecedented angular and 
energy resolutions, has produced the most detailed maps of the gamma ray sky for a wide 
range of energies. Intriguingly, as was reported by several independent groups [1-9] and 
also by Fermi Collaboration itself [10], the Fermi-LAT data have revealed the presence of 
an extended excess of gamma rays over the modeled foreground and background emissions 
towards the Galactic Center (GC). Although several astrophysical mechanisms, such as the 
thousands of unresolved millisecond pulsars [11-13] and the interactions between comic rays 
(CR) and interstellar gases [14-18], have been proposed to interpret this Galactic Center 
Excess (GCE), they usually fail to generate the morphology and energy spectrum of the 
GCE simultaneously 1 . So in this work, we instead consider another possibility that the 
GCE is produced by the annihilation of DM. Although this interpretation has been con¬ 
strained by the measurements of CR such as the Fermi-LAT detection of the gamma-rays 
from dwarf spheroidal galaxies (dSphs) [20-23], the non-observation of spectral features 
in the AMS-02 measurements of CR positron [24-27], and PAMELA observations of the 
CR anti-protons [28-34] , it still remains a most attractive one not only because the excess 
emission shows spectral and morphological properties consistent with a telltale sign from 
DM annihilation, but also because in such an interpretation, the annihilation cross section 
required to explain the GCE is of the right size to account for the DM density from thermal 
freeze-out. 

So far there have been a large number of attempts to explain the GCE by DM anni¬ 
hilation in various new physics models [35-102], In the early analyses of the annihilations, 
great efforts were focused on the channels XX ~I bb with m ^ ~ 35GeV and XX T 7 with 
m% ~ lOGeV since they can reproduce well the GCE spectrum obtained at that time. 
Recently a critical reassessment of the DM interpretation was made by examining in a 
comprehensive way the foreground and background uncertainties [9]. It was found that 
taking the estimated uncertainty in the high-energy tail of the spectrum into account, a 
much larger number of DM annihilations are able to fit well the 7 -ray data than previously 
noted [73, 79]. Explicitly speaking, as far as the annihilation XX ~t bb is concerned, now 
the mass of DM is extended to a broader range from 30GeV to 70GeV in explaining the 
GCE [73, 79]. Other annihilation channels such as DM annihilation into light quark pairs 
and even gluon pair are also able to provide a good fit to the GCE [79]. More strikingly, 
this new analysis opens up a very good solution usually neglected before, namely DM an¬ 
nihilation into a pair of light non-standard Higgs bosons [49, 50] . This important progress 

1 An exception may be the mechanisms recently proposed in [19]. 
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motivates us to renew the solutions to the GCE in supersymmetric theories, which usually 
predict the lightest neutralino x§ as a natural DM candidate. 

As the most economical realization of supersymmetry, the Minimal Supersymmetric 
Standard Model (MSSM) is unsatisfactory in explaining the GCE due to the following four 
reasons [ 66 , 93]. First, the relic density of DM has required its mass to be larger than 
about 40GeV [103]. In this case, the annihilations XiXi —> Tf,qq with q denoting a light 
quark can not provide a good fit any more. Second, except for excessive fine-tuning cases 
the LHC experiments have pushed the lower mass bounds for the CP-odd Higgs boson 
and the bottom squarks up to several hundred GeV. As a result, the cross section of 
DM annihilation into bb in present day is too small to significantly contribute to the GCE 
[ 66 , 93]. Third, due to the small velocity of DM in our galaxy, the annihilation rate for 
DM into SM-like Higgs pairs is p -wave suppressed. Consequently this channel is not large 
enough to generate the GCE. Finally, as for the annihilations XiXi ~I WW. , ZZ , their fits 
to the GCE spectrum indicate that regardless of their annihilation rates the corresponding 
p- values are always less than 0.04 [73, 79]. This means that the annihilations can not 
generate the proper spectrum shape for the GCE. We note that for a given parameter 
point of the MSSM, DM usually annihilates into multiple final states. In this case, the 
situation can not be improved greatly because, due to the particle spectrum of the MSSM 
allowed by the current experiments, either the total cross section falls short for the GCE, 
or the dominant annihilation channel can not reproduce the GCE spectrum well [92]. 

Given the problems of the MSSM, we consider to interpret the GCE in the Next-to- 
Minimal Supersymmetric Standard Model (NMSSM) with a Z 3 symmetry, which is the 
simplest gauge singlet Higgs extension of the MSSM [104] . Distinguished from the MSSM, 
the NMSSM predicts three singlet-dominated particles: one neutralino, one CP-even and 
one CP-odd Higgs bosons. These particles are rather special in that all of them can be 
simultaneously lighter than about 100 GeV, and that the couplings for the interactions 
among themselves are determined by the parameter k, which alone is able to predict the 
right rates for some annihilation channels to explain the GCE (see the following discussion). 
These features make the NMSSM with a singlet-dominated DM well suit to account for 
the GCE because, as we will show below, some golden channels for the GCE need light 
particles to act as the DM, the mediator and/or the annihilation final state. 

We note that the interpretations of the GCE in the NMSSM have been intensively 
discussed in [57, 58, 66 , 68 , 69, 93]. However, in [57, 58, 66 , 68 , 69] the authors did not 
consider the systematic uncertainties mentioned above. As a result, the model parameter 
space they considered is much narrower than that of this work and the obtained conclusions 
were incomplete. While for [93], although the authors have taken the uncertainties into 
account, they considered the parameter space characterized by a large A which is different 
from our discussion. 

The aim of this work is to explore any possible solution to the GCE in the Z 3 NMSSM. 
For this end, we perform a sophisticated scan over the model parameters by considering 
various experimental constraints such as the DM relic density, the Higgs data as well as 
the observation of dwarf galaxies. We use the correlation matrix presented in [9] to include 
the systematic uncertainties on the GCE spectrum and only keep the parameter points 
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that can reproduce well the spectrum. In our study we mainly consider a singlino-like DM 
which is believed to interpret the GCE without excessive fine tuning. As we will show 
below, the annihilation X 1 X 1 HiA\ with HiA\ denoting a scalar-pseudoscalar Higgs pair 
may provide the best fit to GCE, and the canonical annihilation X 1 X 1 ~> bb still remains 
a satisfactory solution except that m^o is now allowed to vary within a broader range. 
Moreover, it is interesting to see that the mixed annihilation into W + W~ and H{A\ final 
states is also able to generate a spectrum consistent with the GCE. These conclusions are 
quite different from previous studies in the NMSSM. 

This paper is organized as follows. In Section II, we introduce some of the characteristic 
features of NMSSM, the basic knowledge about the GCE and our strategy for the parameter 
scan. In Section III, we discuss in detail the interpretations of the GCE when H 2 is the 
SM-like Higgs boson, and in Section IV, we carry out a similar study but for the case that 
Hi acts as the SM-like Higgs boson. We draw our conclusion in Section V and provide 
more information of the NMSSM couplings in the Appendix. 

2 Fitting the GCE in the NMSSM 

2.1 Theoretical setup for the GCE in the NMSSM 

We start our analysis by recapitulating the basics of the NMSSM. As one of the most eco¬ 
nomical extensions of the MSSM, the NMSSM introduces one gauge singlet Higgs superfield 
in its matter content, and since one purpose of the extension is to solve the //-problem of 
the MSSM, a Z 3 symmetry is usually adopted in the construction of the superpotential to 
avoid the appearance of parameters with mass dimension. As a result, the superpotential 
of the NMSSM and the soft breaking terms in Higgs sector are given by [104] 

ITNMSSm = Wf + A H u ■ H d S + 1 kS 3 , (2.1) 

3 

Koft SSM = + rn 2 d \H d \ 2 + m 2 s \S\ 2 + ( \A X SH U ■ H d + i kA k S 3 + h.c.), (2.2) 

where Wp is the superpotential of the MSSM without the //-term, H u , H d and S are Higgs 
superfields with H u , H d and S acting as their scalar components respectively, the dimen¬ 
sionless coefficients A and n parameterize the strengthes of the Higgs self couplings, and 
rh u , in d , ins, A\ and A K are soft-breaking parameters. In practice, after the electroweak 
symmetry breaking the soft-breaking squared masses rh^, in^ and m 2 are traded for mz, 
tan/3 = v u /v d and // = Xv s as theoretical inputs. 

Due to the presence of the superfield S, the NMSSM contains a singlino held which is 
the fermion component of S, and one more complex Higgs held S compared to the MSSM. 
As a result, the neutralino mass eigenstates Xi (with i ranging from 1 to 5) are the mixtures 
of bino, wino, higgsinos and singlino, and the CP-even (odd) Higgs mass eigenstates Hi 
with i = 1,2,3 (Ai with i = 1,2) are mixtures of the real (imaginary) parts of H u , H d 
and S. Throughout this paper, we assume the mass order m^o < m^o < ■ ■ ■ < m^o for 
neutralinos, and to//, < m # 2 < rn h 3 , mA 1 < mA 2 f° r Higgs bosons. 

There are three distinguished features in the NMSSM. One is that DM in the NMSSM 
may be either singlino-dominated or bino-dominated. As expected, the properties of a 
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singlino-dominated DM are quite different from those of a bino-dominated DM, which 
makes the DM physics in the NMSSM much richer than that in the MSSM [105]. Another 
feature is that, in the presence of a singlino-dominated DM with mass below lOOGeV, the 
singlet-dominated CP-even and CP-odd Higgs bosons can be simultaneously lighter than 
about lOOGeV [105, 106], and the strengthes for the interactions among these particles are 
determined by the parameter n which may be as large as 0.1. This feature, as we will show 
below, makes the NMSSM with a singlino-dominated DM well suit to explain the GCE. 
In the appendix, we list the properties of these particles used in our analysis. The other 
feature is that either H\ or H 2 in the NMSSM can act as the SM-like Higgs boson [107] 
and generally speaking, H -2 as the SM-like boson is more attractive from phenomenological 
point of view and also from naturalness argument. 

In the DM explanation of the GCE, the observed 7 -ray originates mainly from the 
cascade decays of the annihilation final states. In the NMSSM, the possible annihilation 
final states include //, VV, H L Hj , A{Aj and HiAj [108], where / (V) denotes any of the 
fermions (vector bosons) in the SM, and Hi ( Aj ) denotes a CP-even (CP-odd) Higgs boson. 
In this work, we are particularly interested in the annihilations X 1 X 1 ~W + W~, Hj A\. 
These annihilations proceed through the s-channel mediator of a Z boson or a Higgs boson 
with an appropriate CP quantum number, and also proceed through the i/u-channel ex¬ 
change of a sbottom, a chargino and a neutralino respectively. The complete expressions 
of the annihilation cross sections are rather complicated, but in non-realistic limit, i.e. the 
velocity of DM approaching zero, some contributions become unimportant. In this case, 
the velocity weighted annihilation cross section can be approximated by [108] 
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where Cxyz denotes the coupling of the interaction involving the particles X, Y and Z, 
is the width of the CP-odd state Aj, oj = fej = m 2 ± /m‘l v , f l j J is 

the coupling coefficient for X? (Xixt R W T ) interaction, and 5 = (2?n^o — ( run, + 
m Al ))/ 2m^o. In getting Eq.(2.5), we note that a good fit to the GCE requires that the 
HiA\ final state is produced close to threshold, i.e. 5 ~ 0, so we can expand (< THiA^o 
in terms of 5. Then the first two terms on the right hand of Eq.(2.5) come from the left 
diagram of Fig.l, and the last term comes from the right diagram of Fig.l. 

The flux per unit solid angle at some photon energy Ey, which is observed by Fermi- 
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Figure 1. Feynman diagrams contributing to the annihilation X 1 X 1 —^ H, A[ with x® (j from 1 to 
5) denoting any of the five neutralinos. A u-channel diagram in associated with the t-channel one 
is assumed. 
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where dN YY /dE is the photon spectrum generated by the annihilation XiXi XY, pdm 
is the DM profile and the integral over p 2 DM is along the light-of-sight (LOS) at an angle 6 
towards GC. In the DM interpretation of the GCE, a generalized Navarro, Frenk & White 
(NFW) DM profile is usually adopted, and its expression is given by [109, 110] 


P(r) = Po 



( f + r e/Rs 
V 1 + r/R s 


3-7 


(2.7) 


with slope parameter 7 = 1.26, scale radius R s = 20 kpc and the local DM density 
Po = 0.4GeV/cm 3 at the radial distance of the sun from the galactic centre Tq. Here 
the coordinate r is centered on the galactic centre and can be expressed as r 2 (s,6) = 
r|, + s 2 — 2rQSCos0 with s and 9 being the LOS distance and the aperture angle between 
the axis connecting the earth with the galactic centre and the LOS respectively. 

In our study, we use the package micrOMEGAs-3.6.9.2 [111] to calculate the DM relic 
density and with the help of PYTHIA [112] to generate the flux in Eq.(2.6). Note that 
in any explicit model, DM usually annihilates into multiple final states. In this case, the 
different fluxes are summed over. 


2.2 Parameter scan strategy for GCE solution 

We simplify our scan over the NMSSM parameter space by fixing the parameters that are 
not closely related to the DM studies. The soft SUSY breaking parameters in the squark 
sector are all fixed to be 2 TeV except that we vary those for the third generation to 
generate a CP-even Higgs near 125 GeV. We assume A t = A \, and Mjj 3 = Md 3 to reduce 
the number of free parameters. Similarly, all of the soft SUSY breaking parameters in 
the slepton sector are fixed to be 300 GeV to explain the discrepancy of the measured 
value for muon anomalous magnetic moment from its SM prediction. As for the gaugino 
sector we abandon the Grand Unified Theory relation and fix the wino mass and gluino 
mass at 1 TeV and 2 TeV respectively. Consequently, the remained free parameters include 
tan/3, p, A, k, A\, A k in the Higgs sector, Mq 3 , Mjj 3 and At for third generation quarks 
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and the bino mass M\, which are all defined at the scale of 2 TeV in the scan. We use 
NMSSMTools-4.3.0 [113] to scan intensively the following NMSSM parameter region: 


1 < tan/3 < 40, 0 < A < 0.7, 0 < |k| < 0.7, |Mi| < 600 GeV, 

0 < A\ < 5 TeV, \A k \ < 2 TeV, \At\ < 5 TeV, 

100 GeV < /j, < 600 GeV, 200 GeV < mQ 3 ,mu 3 < 5 TeV. (2-8) 


Table 1 . Favored parameter region of the NMSSM to explain the GCE, which are classified by the 
dominant final state in DM annihilations. These annihilations are called Solution I, II, III, IV and 
V respectively in the following discussion. All input parameters are defined at 2TeV and quantities 
with mass (annihilation cross section) dimension are in unit of GeV (10 -26 cm 3 /s). 



H 2 is SM-like 

Hi is SM-like 

bb 

AT, M, 

W+W~ 

Hi A! 

H 2 Ai 

tan(3 

(10,30) 

(8,40) 

(8,40) 

(15,20) 

(12,20) 

Ml 

(-600,-60) 

(-600,-80) 

(-600,-130) 

(-130,-90) 

(-200,-100) 

A 

(0.2, 0.7) 

(0.2, 0.6) 

(0.2, 0.4) 

(0.6, 0.7) 

(0.4,0.7) 

Av 

(0.02,0.12) 

(0.07,0.15) 

(0.09,0.14) 

(0.10,0.14) 

(0.11,0.16) 

4 

(160,300) 

(110,210) 

(110,160) 

(220,270) 

(210,270) 

A\ 

(2400, 5000) 

(830, 5000) 

(970, 5000) 

(3900, 5000) 

(2900,5000) 

a k 

(-210,-70) 

(-60,22) 

(—70,10) 

(-65,-16) 

(-66,5) 

At,b 

(-4300, 3900) 

(-4600,4700) 

(-4700, 3900) 

(-2200,2000) 

(-3400,4000) 

m q 3 

(300, 5000) 

(350, 5000) 

(500, 5000) 

(1200,4600) 

(700,4800) 

Mu 3 ,d 3 

(250, 5000) 

(270, 5000) 

(400, 5000) 

(250,5000) 

(1400,5000) 

m Hl 

(15,102) 

(61,119) 

(83,110) 

(124,127.4) 

(124.5,127.4) 

m,H 2 

(122.8,127.8) 

(122.7,128) 

(123,128) 

(125.6,142) 

(125.7,146) 

(av ) 0 

(0.17,1.9) 

(0.29,1.8) 

(0.44,1.6) 

(0.34,1.2) 

(0.38,1.5) 

m^o 

Xi 

(31,70) 

(62,114) 

(84,102) 

(71,87) 

(80,127) 

m-o 

(67, 298) 

(83, 233) 

(114,165) 

(86,128) 

(100,192) 

m x± 

(166,297) 

(117,214) 

(117,158) 

(225, 270) 

(218,266) 

m Al 

(58,133) 

(9,109) 

(10,105) 

(9,38) 

(16,95) 

m H ± 

(3510,4666) 

(1477, 3447) 

(2818,2968) 

(3883, 4941) 

(2928, 4740) 

Xmin 

(23,35) 

(21,35) 

(24,35) 

(21,35) 

(21,35) 


The process to retain the parameter points include the following steps: 

• We require the DM to be singlino-dominated and satisfy m^o <150 GeV, and impose 
all the experimental constraints encoded in NMSSMTools-4.3.0 [113] which include 
the relic abundance at 3cr level (0.107 < H/i 2 < 0.131), LUX exclusion bound at 90% 
C.L., various B-physics measurements as well as the discrepancy of muon magnetic 
moment at 2cr level. We also consider various electroweak precision data calculated 
in [114]. 

• We consider the constraints on the Higgs sector with the package HiggsBounds-4.1.2 
[115] which contains the data from LEP, Tevetron and LHC. For the SM-like Higgs 
boson, we further perform a fit to the data with the package HiggsSignal [116] and 
keep the 2 er samples. 
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• We use micrOMEGAs-3.6.9.2 [111] to calculate the DM annihilation cross section at 
present day, and then impose the constraints from dSphs by the data in [23] for the 
bb annihilation channel and with the method introduced in [93] for the HiA\ final 
states. 

• We also use micrOMEGAs-3.6.9.2 [111] to generate the 7 -ray spectrum. Considering 
the astrophysical uncertainties which may come from the errors in our setting on the 
local DM density pq, the scale radius R s and the inner slope parameter 7 in Eq.(7), 
for each parameter point we allow an uncertainty factor A in the range of (0.17, 5.3) 
for the annihilation cross section, or equivalently for the height of the gamma-ray 
spectrum in Eq.(2.6) [79]. Then for the A —tuned 7 -ray spectrum, we perform a 
fit to the residual GCE spectrum obtained in [9] by using the publicly available 
covariance matrix, which include both the statistical and systematic uncertainties of 
the measured flux. The corresponding Xsp function is calculated by [9, 79]: 



where is the covariance matrix, dN/dEi is the measured flux in the f-tli energy bin, 
and dN/dEi is the flux predicted by the NMSSM, which depends on the parameter 
point and also on the factor A. 

We define the GCE y 2 as the minimum value of Xsp(-A) among different choices of 
A, Xgce = m i n (Xsp(-A))i an d keep the parameter points that satisfy Xgce — 35.2. 
These points are assumed to have the capability to explain the GCE at 95% confidence 
level for 23 degree of freedom [9]. 

The parameter ranges of the GCE solutions are listed in Table.l, which are classified 
by the dominant final state in DM annihilations (see the following discussion). For the first 
three types of the DM annihilations H 2 acts as the SM-like Higgs boson, while for the last 
two types H\ corresponds to the SM-like Higgs boson. One distinguished feature that Table 
1 exhibits is that all the singlet dominated particles in the GCE solutions, including DM, 
the singlet-dominated CP-even and CP-odd Higgs bosons, are lighter than about 150GeV. 
This feature, as we will emphasized below, makes the NMSSM well suit for explaining the 
GCE. 

3 GCE solutions with H 2 being the SM-like Higgs boson 

In this section, we exhibit the features of the GCE solutions for the case that DM is 
singlino-dominated and H 2 acts as the SM-like Higgs boson. All the solutions considered 
in this work survive the constraints listed in last section and meanwhile can explain the 
GCE at 95% C.L.. 

In Fig.l we project the solutions on (av)o — m^o plane (upper panel) and Xgce ~ m \\ 
plane (lower panel). Solutions marked by red square, blue triangle and black asterisk 
correspond to the cases that DM annihilates with the largest branching ratio into bb, H\A\ 
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Table 2. Detailed information of the benchmark points used in our discussion. Quantities with 
mass, annihilation and scattering cross section dimension are in unit of GeV, cm 3 /s and pb respec¬ 
tively. 


Point 

tan/3 

A 

K, 



a k 

Ad 3 ,u 3 

Mi 

Mq 3 

Mu 3t D 3 

PI 

16 

0.36 

0.04 

241 

3891 

-136 

420 

-472 

4127 

4445 

P2 

12 

0.46 

0.12 

179 

2036 

-6 

-2354 

-209 

2197 

3673 

P3 

13 

0.27 

0.11 

130 

1899 

-5 

-524 

-170 

4098 

4384 

P4 

18 

0.69 

0.12 

243 

4518 

-43 

-320 

-103 

1436 

4308 

P5 

17 

0.66 

0.13 

226 

3923 

-17 

1138 

-97 

4540 

1286 

P6 

18 

0.66 

0.15 

217 

4048 

-24 

2050 

-103 

4170 

1452 

P7 

15 

0.50 

0.13 

255 

4085 

-35 

2621 

-131 

2935 

4468 

Point 

m Hl 

m H2 

m H ± 

m,A l 

m^o 

m^o 

m~± 

Xi 

Sr (h2->x°x?) 

Br( h 1-1A1A1) 

S r (ft2^AiAi) 

PI 

40 

125 

3960 

99 

50 

256 

248 

0.54% 

0 

0 

P2 

99 

125 

2065 

66 

87 

178 

183 

0 

0 

0 

P3 

99 

126 

1823 

48 

92 

126 

134 

0 

88.51% 

7.06% 

P4 

126 

133 

4452 

20 

78 

102 

249 

0 

6.50% 

95.32% 

P5 

125 

126 

3883 

27 

81 

96 

231 

0 

4.86% 

94.01% 

P6 

126 

129 

4022 

33 

85 

101 

222 

0 

4.55% 

95.74% 

P7 

125 

145 

4068 

69 

121 

127 

262 

0 

0 

94.11% 

Point 

2 

Xgce 

p-value 

(av) 0 

{<jv)\ Tf 

SI 

(Tp 

_SD 

(Tp 

^bb 

Rfj 2 Ai 

1 

+ 

RH2A1 

PI 

23.3 

0.44 

6.1E-27 

2.8E-26 

2.3E-15 

1.5E-04 

90.4% 

0.0% 

0.0% 

0.0% 

P2 

21.6 

0.54 

1.4E-26 

2.6E-26 

8.2E-10 

1.3E-03 

0.0% 

96.8% 

3.1% 

0.0% 

P3 

24.8 

0.36 

1.3E-26 

2.7E-26 

1.9E-10 

8.5E-04 

0.1% 

46.1% 

47.2% 

0.0% 

P4 

22.4 

0.50 

8.6E-27 

2.8E-26 

5.5E-10 

1.6E-03 

0.1% 

94.3% 

0.1% 

5.4% 

P5 

21.4 

0.55 

9.9E-27 

3.1E-26 

4.4E-10 

1.7E-03 

0.1% 

68.3% 

0.2% 

31.4% 

P6 

21.6 

0.54 

8.3E-27 

3.2E-26 

9.8E-10 

2.0E-03 

0.1% 

42.1% 

2.1% 

55.5% 

P7 

23.7 

0.42 

7.9E-27 

3.3E-26 

1.1E-09 

5.5E-04 

0.1% 

2.2% 

7.2% 

84.0% 


and W + W~ final states respectively, which hereafter are collectively called Solution I, 
Solution II and Solution III correspondingly. Then the upper panel of Figure. 2 indicates 
that, for the ranges 30GeV < m^o < 40GeV, 50GeV < m^o < 62GeV and 63GeV < 
m^a < 70GeV, Solution I is viable, while for 63GeV < m^o < 115GeV and 83GeV < 
m^o < lOOGeV, Solution II and Solution III can account for the GCE respectively. For any 
of the solutions, the (crv)o is larger than 1.7 X 10 -27 cm 3 /s, and its lower bound increases 
monotonically as Xi becomes heavier. The reason for the latter behavior is that, for a 
heavier DM, its number density is smaller. So to obtain the same photon flux for the GCE, 
a larger cross section is needed. 

The lower panel of Fig. 2 indicates that the best interpretation in Solution I comes 
from m^o ~ 50GeV with x%ce — 23 and a p-value of 0.44. This conclusion coincides 
with that of [79], which was obtained in a model independent way and for a pure bb 
annihilation channel. For Solutions II and Solutions III, the best interpretations locate at 
m^o — 87GeV with Xgce — 21-6 and m^o — 92GeV with Xgce — 24.7 respectively, and 
the corresponding p-values are 0.54 and 0.36. These two solutions, within our knowledge, 
were rarely discussed in previous literatures about the NMSSM. Moreover, we checked 
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m- x o{GeV) 

Figure 2. The GCE solutions of the NMSSM for a singlino-dominated DM with H 2 acting as the 
SM-like Higgs boson, which are projected on the (av) o — m^o plane (upper panel) and Xgce ~ 
iTiyO plane (lower panel). Solutions marked by the red square, the blue triangle and the black 
asterisk correspond to the case that DM annihilates in present day mainly by the channels XiXi _■ y 
bb, H\Ai,W + W~ respectively, which are collectively called Solution I, II and III correspondingly. 

that, in the case of m^o ~ rnz /2 (m^.o — 62GeV), DM annihilated in early universe mainly 
through a nearly on-shell Z boson (SM-like Higgs boson). Since nowadays this dominant 
annihilation channel is helicity (p- wave) suppressed, (av) o can not reach the size required 
for the GCE. 

In Table 2, we present detailed information of three benchmark points PI, P2 and 
P3 for Solution I, II and III respectively. This table indicates that the sole annihilation 
channel XiXi ~t bb or XVxV H\A\ can be responsible for the GCE; while for the channel 
XiX? —> W + W , it must mix sizeably with the channel XiXi H\A\ to account for the 
GCE. We will return to this issue later. 

In our calculation, we found that the condition on the GCE x 2 can reduce the number 
of the parameter points that survive the constraints by more than 90%. This implies 
that the GCE has non-trivial requirements on the parameters of the NMSSM, especially it 
suggests that some of the independent parameters may be correlated. Motivated by this 
thought, we study the correlations among the parameters A, k, //, m^o and m^ 1 which 
are important parameters in the interpretation of the GCE and show the corresponding 
results in Fig. 3. In the following, we concentrate separately on each kind of the solutions 
and investigate its features. Such a study is helpful to understand the correlations in Fig.3 
and also the properties of the benchmark points listed in Table. 2. 
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3.1 Solution I - the bb annihilation channel 


Among the solutions to the GCE, Solution I is the most intensively studied one. After 
considering the systematic uncertainties, one important improvement of Solution I over its 
previous version is that DM mass is now allowed in the range from 30GeV to 70GeV, which 
is much wider than before. 

The key features of Solution I are as follows: 


• The lighter CP-odd Higgs boson is correlated with DM by mA 1 — 2m^o. This corre- 

1 A-l 

lation is shown in the upper right panel of Fig. 3 which means that the annihilation 
proceeds resonantly. 


This feature can be understood as follows. In Solution I, the heavy CP-odd Higgs 
boson is doublet-dominated with its mass usually at TeV scale. Then Eq.(2.3) in¬ 
dicates that the main contribution to the annihilation comes from the moderately 
light A\, which is singlet-dominated. With the formula presented in Eq.(A.20) and 
v s = n/\ > 450GeV shown in the lower left panel of Fig.3, one can get 


c Mm c U - * 2 * 2 < - 


mb V 

» ) 


< 


(AY 

\450 ) 


(3.1) 
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This inequation means that the couplings involved in the annihilation are highly sup¬ 
pressed so that the process must proceed resonantly to ensure (o bb v)o 10 26 cm 3 /s. 
Moreover, our results indicate that the width of A\ is very small, T^ < 10“ 2 MeV. 
So as rriA 1 approaches 2m^o, the denominator in Eq.(2.3) tends to vanish and a small 
k in Eq.(3.1) is then suffice to predict the right rate of the annihilation for the GCE. 
This character is illustrated in the upper left panel of Fig. 3. In fact, a small k is also 
favored to predict light Xa and A]. which can be seen from Eq.(A.3) and Eq.(A. 8 ). 

• The parameter [i is upper bounded by about 300GeV, which is shown in the lower 
panels of Fig. 3. 

This feature is actually required by the DM relic density [69]. Generally speaking, in 
order to predict the measured H/? 2 , the velocity weighted cross section (crv) should 
be around the canonical value 3 x ICE 26 cm 3 /s at freezing out (see for example 
points in Table 2). Since 2m. x o/mA 1 > 1 in Solution I, (crv) for the annihilation 
XiXi —> A\ —> bb at present day is usually larger than that at freezing out due to 
the thermal broadening [117]. Since the dwarf galaxy measurements have required 
( cr bb v)o <2x 10 -26 cm 3 /s (see Fig.2), new contributions such as those mediated by a 
Z boson or a CP-even Higgs boson must intervene for the DM annihilation in early 
Universe, and a moderately small /r can accelerate the annihilation [69]. 


• Solution I suffers from severe fine tuning problem. Explicitly speaking, beside the cor¬ 
relation rn a \ — 2 m^o, there exits another strong correlation observed in our analysis, 
which is given by 


m x o/GeV ~ 


51 - 475k, 
37 + 325k, 
49 + 175k, 


for 30GeV < m^o < 40GeV or 0.024 < k < 0.045, 
for 50GeV < m^o < 62GeV or 0.038 < k < 0.07, 
for 63GeV < m+o < 70GeV or 0.08 < k < 0 . 12 . 

A-l 


These correlations make Solution I in the NMSSM quite unnatural to explain the 
GCE. 


• We checked that Br(A\ -x 77 ) < 5 X 10 4 so that the 7 -ray spectral line generated 
by XaXa —> A\ —> 77 is suppressed. 

• Since ^ 60GeV for most cases in Solution I, the SM-like Higgs boson H -2 may 
decay into Xi pair. We checked that Br(H -2 -X X 1 X 1 ) < 18%, which is required by 
the Higgs data at the LHC. 


3.2 Solution II - the H\A\ annihilation channel 

Solution II is quite similar to the interpretations presented in [49, 50, 83, 93, 99, 100], 
which utilize the process XX ftifa ~> / 1 / 1 / 2/2 (4> 1 and ^2 denote scalar or pseudoscalar 
particles, and f\ and fi are SM fermions) for the GCE. These interpretations, as were 
emphasized by the proposers, can easily escape the constraints from DM detection experi¬ 
ments and have been paid more and more attention recently. 

The features of Solution II are as follows: 
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• The singlet-dominated particles satisfy 60GeV < m^o < 115GeV, lOGeV < mA 1 < 
HOGeV, 60GeV < mn 1 < 120GeV and <5 < 0.2, and for most samples there exist 
following relations mn 1 ^ m x° 1 ^ m A 1 - Given n ~ 0.1 which is required to predict 
the right size of the annihilation XiXi HiA\ for the GCE (see below), the particle 
spectrum limits parameters such as A, g and A K in certain regions (see the expressions 
of the tree level masses in Appendix), which are given in Table 1, and also shown in 
Fig-3. 

Note that g is below about 200GeV. In this case, the higgsino-dominated neutralinos 
Xi may decay dominantly into Xi^i instead of into X\Z since the kinematics is 
forbidden. In this case, the LHC search for electroweakinos by trilepton +Eip' lss 
signal is less efficient in ruling out the light higgsinos 2 . Also note that the parameters 
A and g are related by g/GeV « 60 + 260 A for A varying from 0.2 to 0.6 (see lower 
left panel of Fig. 2), which means that v s = g/X > 360GeV. This ensures that the 
expansions for the masses and couplings in Appendix by the power of Xv/g are good 
approximations. 

• The s-channel contributions to the annihilation rate {(Th 1 Ai'u)o in Eq.(2.5) are usually 
much smaller than those from the t/u channel, and among the t/u channel contribu¬ 
tions, the one induced by the exchange of Xi is f ar dominant. As for the contributions 
induced by the two higgsino-like neutralinos, each of them may be sizable, but since 
they cancel each other, the net higgsino contribution is not important. These char¬ 
acters can be understood by the following approximations (see Eq.(A.20)) 


for Higgsino — like Xi and m x° < 0, 
for Higgsino — like and m^o > 0, 

and by the fact that k ~ 0.1 is enough to predict the y,i contributed {(Th.a^o at 
the order of 1CP 26 err? / s (see equation (3.20) in [93]). 

• Since tha 1 < 60GeV for most cases in Solution II (see upper right panel of Fig.2), 
the SM-like Higgs boson H 2 may decay into A\A\ with a sizeable fraction. Given 
that A\ decays dominantly into bb, this will result in 4 b signal for the SM-like Higgs 
boson. We checked that Br(H 2 —> A\Ai) < 24%, where the upper bound comes from 
the constraints of the LHC Higgs data. 

• Since a good fit to the GCE requires that H\A\ is produced close to threshold, the 

annihilation ~> H\A\ will produce spectral line or box-shaped spectrum in 

7 -ray [83, 93]. We checked that Br(H\ -A 77 ) < 1 x 10 -3 for most samples and 

2 In doing [105], we once confronted with the situation quite similar to what we are facing now. Our 
detailed simulation at that time indicated that the trilepton constraint on SUSY is very weak. Moreover, 
in comparison with the case discussed in [118], we find that our case is more difficult to detect since the 
signal is smaller. 


C A 1 x 0 1 X° 1 C H 1 x 0 1 X°i ~ 2 * k2(1 + 2 ~) 2 

r 


C A 1 xlx 0 i C H 1 x° 1 x° i - 


i \ 2 v 2 
4%W 


sin 2 /3, 
sin 2 /?, 
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Figure 4. The GCE y 2 as a function of DM mass for Solution III where DM annihilates mainly 
into WW final state. Here RhxA-l denotes the branching ratio of the annihilation into H\A\ final 
state. 


Br{A\ —> 77 ) < 4 x 1CP 4 for all samples. So current results of the Fermi-LAT search 
for spectral lines [119] can not impose tight limit on Solution II (see [83] for a detailed 
discussion). 

3.3 Solution III - the W + W~ annihilation channel 

In general, the pure annihilation Xx'Xa ~> W + W~ is unable to explain the GCE quite well 
[73, 79], but if it mixes sizably with other annihilation channels, the generated spectrum 
may be improved significantly to account for the GCE. Solution III in the NMSSM belongs 
to this case. 

Solution III has the following features: 

• The W pair must be produced close to threshold to account for the GCE, which 
means 85GeV < m^o < lOOGeV (see right panels of Fig.3). 

• From the expression of (crwwv )0 in Eq.(2.4), one can learn that, if the wino is de¬ 
coupled, the annihilation rate is determined by the higgsino-dominated chargino. In 
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this case, we have 


, 9 AT 9 Au . r- . nv 

fi,L - - 7=^14 ^ -7= —sm/3 ~ v 2g sin (3 -, 

V 2 V 2 M m^o 

/i,R - —~ —\/ 2(7 sin/3—. (3.2) 

v 2 /4 

In getting these expressions, we note v s = ^/X > 400GeV (see lower left panel of 
Fig. 2), and expand IV 13 and IV 14 in terms of Xv//i (see Appendix). We also use the 
approximation 0 ~ 2k///A. Then (crwwv )0 ~ 10~ 26 cm 3 /s and 0 ~ 90GeV 
limit tightly the ranges of the parameters A, k and //, which are shown in Table 1 
and Fig.3. 

Note in Solution III, the parameter /q or equivalently the masses for the higgsino- 
dominated chargino and neutralinos, is less than about 150GeV. Since the splitting 
between a and ra { 0 is less than about 50GeV, such a low value of u is still allowed 

A-l 

by the LHC search for SUSY (see footnote 2 in our discussion on Solution II). 

• The upper left panel of Fig.3 indicates that the parameters A and n are correlated by 

n ~ 0.03 + 0.3A, for 0.2 < A < 0.4. (3.3) 

As a result, we have m^o — 2/r/3. 

• As we emphasized before, the annihilation X 1 X 1 W + W~ must mix sizably with the 
annihilation X 1 X 1 ~t H\A\ to explain the GCE. This, in return, requires appropriate 
masses for H\ and A\ to improve the 7 -ray spectrum generated by the WW state. 
In Fig. 4, we plot the GCE y 2 as a function of DM mass in Solution III with different 
colors denoting the branching ratio of the DM annihilation into H\A\. This figure 
indicates that, with the increase of the branching ratio, the GCE y 2 tends to decrease. 

4 GCE solutions with Hi being the SM-like Higgs boson 

In this section, we investigate the GCE solutions for the case that DM is singlino-dominated, 
and meanwhile H\ acts as the SM-like Higgs boson. We carry out our study in a way similar 
to what we did in Section III. 

In Fig.5 we project the solutions on (<ru)o — m^o plane (upper panel) and Xgce ~ m x\ 
plane (lower panel). For solutions marked by green lozenge, DM annihilates with the largest 
branching ratio into H\A\. while for those marked by red pentastar, DM annihilates mainly 
into H 2 A 1 . In the following, we call these two kinds of solutions Solution IV and Solution V 
respectively. Fig.5 then indicates that, for 70GeV < 0 < 87GeV, Solution IV can explain 

the GCE quite well with the best explanation coming from m^o ~ 81GeV with Xgce — 21.4 
(corresponding to a p-value of 0.55), and for 80GeV < m^o < 130GeV, Solution V is good 
in accounting for the GCE with the best explanation locating at rn-o ~ 85GeV with 

A-l 

Xgce — 21.6 and a p-value of 0.54. 

Compared with the case that H 2 acts as the SM-like Higgs boson, we find that it is 
more difficult to get the GCE solutions if H 1 corresponds to the SM-like Higgs boson. One 
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m^o(GeV) 

Figure 5. Similar to Fig. 2, but showing the solutions for the case that Hi acts as the SM-like 
Higgs boson. For these solutions, DM may mainly annihilate into H\A\ final state (called Solution 
IV in our discussion) or into H 2 A\ final state (Solution V). 


important reason is that the spectrum of the singlet-dominated particles for Solution IV 
and V has non-trivial requirements on the NMSSM parameters, which can not be easily 
satisfied due to the structure of the NMSSM itself. A good example about this argument is 
that we do not find any solutions where DM mainly annihilates into bb. This is due to the 
fact that, given a singlino-dominated DM with 30GeV < m^~o < 70GeV and meanwhile a 

A-l 

singlet-dominated A\ satisfying — 2m^o, the singlet-dominated CP-even Higgs boson 
is usually lighter than the SM-like Higgs boson [69]. 

In Table 2, we present detailed information for benchmark points P4, P5, P6 and 
P7 with points P4 and P5 belonging to Solution IV and points P6 and P7 belonging to 
Solution V. This table shows that for 80GeV < m^o < 86GeV, DM may annihilate into 
H\A\ and H 2 A\ states with comparable rates to explain the GCE (see points P5 and P6), 
while for m^o ~ 78GeV (m^o ~ 120GeV), the sole annihilation channel XiXi H\A\ 
(X?X? H 2 A 1 ) can be responsible for the GCE, see point P4 (P7). 

In Fig. 6, we show the correlations among the parameters A, k. /i, m^o and m J 4 1 . 
This figure is supplement to Table 1, and as we will show below, it is helpful for our 
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Figure 6. Similar to Fig. 3, but showing the correlations for Solution IV and Solution V. 


understanding on Solution IV and V. 

4.1 Solution IV - the H\A\ annihilation channel 

Solution IV has the following features: 

• The H\A\ state must be produced close to threshold to explain the GCE, which is 
reflected by 5 < 0.1 from our results. 

• The favored spectrum for the singlet-dominated particles is 71GeV < m^o < 87GeV, 
lOGeV < niA 1 < 40GeV and 126GeV < mn 2 ^5 142GeV. Given n ~ 0.12 which is 
required to predict the right size of the annihilation XiXi H\A\ for the GCE (see 
below), this spectrum limits parameters such as A, k. /r and A K within rather narrow 
ranges, which are given in Table 1 and also shown in Fig. 6. 

Compared with Solution II, we find in Solution IV that, in order to predict a heavier 
singlet-dominated CP-even Higgs boson, the parameter fi usually takes a larger value, 
220GeV < n < 270GeV. As a result, A must exceed about 0.6, which can be inferred 
from the relation m-o ~ 2nfj,/X ~ 2 x 0.12 x ///A ~ 80GeV. This relation also 
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suggests that v s = ft/ A > 300GeV or Xv/g < 0.6, which makes the expansions listed 
in Appendix feasible. 


• Similar to Solution II, the s-channel contributions to the annihilation rate {(Th 1 a 1 'v)q 
in Eq.(2.5) are usually much smaller than those from the t/u channel, and among 
the t/u channel contributions, the one induced by the exchange of x 1 is dominant. 
However, since Hi now is the SM-like Higgs boson (instead of a singlet-dominated 
particle in Solution II), there still exists a slight difference between the two solutions. 
Explicitly speaking, we find that each higgsino contribution to the annihilation is 
comparable in magnitude with the Xi contribution, but since the two higgsino con¬ 
tributions cancel each other, the total higgsino contribution is small. These features 
can be explained by the following formula (see Eq.(A.20)) 


Xv 

C A ,v 0 v°C h ~o~o ~ —4 in 2 —sin 2 d ~ — 8k d — sin 2 /?, 

^lXiXi "1X1X1 I, r m ~n 




m^o 

A1 


(4.1) 


C A 1 x 0 1 X° i C H 1 x < iXi 


\X 2 ^f sin 2 /3, for Higgsino — like xf/ and m^o < 0, 
— ^A 2 ^ sin 2 /3, for Higgsino — like and rri^o > 0, 


and also by comparing Eq.(4.1) with equation (3.20) in [93] to conclude that k ~ 0.12 
is enough to predict the contributed {cth 1 a 1 v)o at the order of 10 ~ 26 cm 3 /s. 


• Since < 40GeV for all cases in Solution IV (see upper right panel of Fig. 6 ), the 
SM-like Higgs boson H\ will decay into A\A\. We checked that Br(H\ -» A\A\) < 
24% as required by the Higgs data. 


• We also checked that A\ -A bb is the dominant decay mode of A\, and Br{A\ —> 
77 ) < 5 x 10 ~ 5 for all samples. 


4.2 Solution V - the HqA\ annihilation channel 

Since H 2 in Solution V is singlet dominated, the features of Solution V should be similar 
to those of Solution II. The differences mainly come from the following aspects: 


• The spectrum of the singlet dominated particles. In Solution V, the favored spectrum 
is 80GeV < m^o < 130GeV, 18GeV < < lOOGeV and 125GeV < mn 2 ^ 

146GeV with < rn^o < mn 2 an d ^ < 0.1. Corresponding to such a spectrum, 
the parameter space of Solution V differs greatly from that of Solution II, which can 
be seen from Table 1 and also from Fig. 6 . 


• The phenomenology of some relevant particles. For example, in both Solution IV and 
Solution V, the favored value of /1 is uplifted in comparison with that in Solution II. 
As a result, the higgsino-dominated neutralinos may decay into Zx ?, which makes 
them to be potentially detected at 14-TeV LHC by trilepton +Etp iss signals [69]. 


- 18 - 




Figure 7. Spin-independent (SI) and Spin-dependent (SD) cross sections for DM-nucleon scattering 
as a function of DM mass. Solutions in this figure are taken from Fig. 2 and Fig. 5 with the same 
symbolic conventions. 


5 Explore the GCE solutions in future DM experiments 

In this section we investigate to what extent the GCE solutions will be explored in future 
DM direct detection experiments such as XENON-1T and LUX experiments [120], which 
will improve current experimental sensitivities to DM-nucleon scattering cross sections by 
up to three orders. In Fig.7, we project our solutions on m^o — ap 1 and m^o — Op D 
planes with a^ 1 and a^ D denoting the spin-independent (SI) and spin-dependent (SD) 
cross sections respectively. The left panels in the figure are the results for the case that 
H 2 acts as the SM-like Higgs boson, and the right panels are those for the case that H\ 
corresponds to the SM-like Higgs boson. The dotted lines, solid lines, dashed lines and 
dash dotted lines are the sensitivities to the cross sections set by the XENON-100, LUX, 
XENON-IT and LZ experiments respectively. Note that so far the XENON-100 experiment 
has imposed constraints on both SI and SD cross sections, while the LUX experiment only 
obtained limits on the SI cross section. 

For c 7 $ 1 in the H -2 case, we can see from Fig.7 that the future XENON-1T experiment 
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is able to probe a large portion of the GCE solutions, and the LZ experiment can test even 
more solutions. Anyhow, there still exist some solutions remaining untouched by these 
future experiments. This conclusion can be understood as follows. In the NMSSM after 
considering the current experimental constraints on sfermion masses, the main contribution 
to (ip J comes from the t-channel process mediated by the CP-even Higgses H\ 2 - In this 
case, the Wilson coefficient f qi for the operator XiXiQiQi is given by [121] 


C H 1 x 0 1 X° 1 C Hiq i q i C H 2 X° 1 x\ C H 2 qiqi 


iqt 


2 rn 2 Hi 


+ 


2m H 2 


(5.1) 


where C Hl ^o^o ~ —\/2k( 1 + 2Au//r) and C Hr> ^ o^.o ~ 2^/2n\v/fi, which are given by 
Eq.(A.20) and Eq.(A.21) respectively. Then Eq.(5.1) indicates that, if k is small or if 
there exists a strong cancelation between the two terms, f qi or equivalently the SI cross 
section will be suppressed. We numerically checked that the untouched solutions has either 
of the two characteristics. 

On the other hand, the story for a^ D in the H 2 case is quite different. From the 
lower left panel of Fig. 7 we can see that the future XENON-1T experiment can test almost 
all of the GCE solutions, let alone the more sensitive LZ experiment. The underlying 
reason is that in the NMSSM with heavy sfermions, the SD cross section gets contribution 
mainly from the t-channel Z-mediated diagram. As a result, the size of the cross section 
is determined by the ZxiXi coupling, which is given by 


9z 


x?x? 


mz 

y/2v 


(7Vi3 - N? 4 ) ~ - 


m z 


AV 


V2v 9 2 


M _ 4k % 

A 2) ~ V2v 


m z X 2 v 2 


P- 


( 1 - 


0 

m xi ■ 


(5.2) 


In getting this expression, we have used the approximations for IV 13 , -/V 14 and m^o. Then 

A1 

from the results presented in Fig. 3, one can infer that except for some rare cases of Solution 
I, the SD cross section is not suppressed too much. 

In a similar way, one can analyze the results for the H\ case, which are shown on the 
right panels of Fig. 7. For example, the upper right panel indicates that the SI cross sections 
in Solution IV and Solution V are usually larger than 10~ 10 p&. This may be understood by 
a weak cancelation between the two terms in Eq.(5.1). Compared with the H 2 case, both 
the SI cross section and the SD cross section in the H\ case are large and consequently, all 
the solutions will be tested by XENON-IT experiment. 

In principle, the GCE solutions in the NMSSM may also be tested by electroweakino 
production processes at the LHC [69]. We will discuss such an issue in our forthcoming 
work. 


6 Summary 

In this work, we took into account the recently reported foreground and background un¬ 
certainties for the GCE and investigated its explanation by DM annihilation in the frame¬ 
work of the NMSSM. We carried out a sophisticated scan over the NMSSM parameter 
space by considering various experimental constraints such as the Higgs data, B —physics 
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observables, DM relic density, LUX experiment and the dSphs constraints. Then for each 
surviving parameter point we performed a fit to the GCE spectrum by using the corre¬ 
lation matrix that incorporated both the statistical and systematic uncertainties of the 
measured excess. Our results indicate that due to the introduction of the gauge singlet 
Higgs superfield, the NMSSM with a singlino-dominated DM has multiple DM annihilation 
channels that are able to explain the GCE quite well, and all of these explanations require 
the singlet-dominated particles (including one neutralino, one CP-even and one CP-odd 
Higgs bosons) to be moderately light. We also discussed to what extent the future DM 
direct detection experiments can explore the GCE solutions, and we conclude that the 
XENON-IT experiment is very promising in testing nearly all the solutions. 

When choosing the scenario of particle spectrum, we focused on a singlino-dominated 
DM and considered the cases that either H -2 or H\ acts as the SM-like Higgs boson. For 
the popular situation that corresponds to the SM-like Higgs, we have the following 
observations on the GCE solutions: 

• The pure DM annihilation channel X 1 X 1 bb or X 1 X 1 can provide a good 

fit to the GCE spectrum, while the channel XiXi W + W~ must mix sizeably with 
the channel XiXi H\A\ to account for the GCE. 

• For the annihilation XiXi bb, DM mass is now allowed in the range from 30GeV 
to 70GeV which is much wider than before. With the help of an appropriate s- 
channel resonance, the singlet trilinear self-coupling parameter k can be as low as 
0.02 to explain the GCE. Moreover, the higgsino mass parameter /j is upper bounded 
by about 300GeV to ensure a correct DM relic density. Since there exist strong 
correlations between independent parameters, such an explanation suffers from a fine 
tuning problem, which is usually less than 1%. 

• The annihilation XiXi H\A\ may provide a better explanation than the channel 
XiXi bb when H\A\ is produced close to threshold, and the best interpretation 
corresponds to a p-value of 0.55. In this kind of explanation, the singlet-dominated 
particles must satisfy GOGeV < o < 115GeV, lOGeV < m a , < HOGeV, GOGeV < 
niH 1 120GeV and 5 < 0.2. This imposes non-trivial constraints on the NMSSM 
parameters, especially that /j must be less than about 200 GeV. Among various 
contributions to the annihilation, the dominant one comes from the y^-contributed 
t/u channel diagrams, in which the parameter k plays an important role in deciding 
the annihilation rate. 

• Apart from the necessary mixing with the H\A\ final states, W + W~ pair in the 
annihilation XiXi W + W~ must be produced close to threshold to account for the 
GCE. A small fi less than about 150 GeV is necessary to increase the annihilation rate 
through the t/u -channel contributions induced by a higgsino-dominated chargino. 
The LHC search for trilepton -\-E™ iSS signal can not exclude such a possibility since 
the electroweakino production rates at the LHC are relatively low, and meanwhile 
since the splitting between u and ?n,~o is compressed. 

A-l 
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• The detection of spin-independent scattering in the future XEN0N-1T and LUX 
experiments are able to cover a large portion of the GCE-favored parameter space, 
while the spin-dependent detection have a stronger potential to test nearly all of the 
relevant parameter region. 

As for the case that H\ acts as the SM-like Higgs boson, the features of the GCE 
solutions are quite different, which are as follows: 

• In comparison with the H 2 case, it is difficult to find GCE solutions when H\ corre¬ 
sponds to the SM-like Higgs boson, and especially we did not find any solution that 
DM annihilates mainly into bb. The reason is, assuming H\ to be the SM-like Higgs 
boson, there must exist sizeable mass splittings among the light singlet-dominated 
particles to explain the GCE, which is difficult to realize in the NMSSM due to the 
theoretical structure itself. 

• For 80GeV < 0 < 86GeV, DM may annihilate into H\A\ and H 2 A 1 states with 

comparable rates to explain the GCE, while for 0 < 80GeV (m^o > lOOGeV), 
the sole annihilation X 1 X 1 H\A\ (xiXi H 2 A 1 ) can be responsible for the GCE. 
For all these solutions, the singlet-dominated particle H 2 and the parameter /x must 
satisfy 125GeV < mn 2 < 145GeV and 210GeV < n < 270GeV. 

• Both the spin-independent and spin-dependent detection in the future XENON-1T 
experiment have a great potential to test the relevant parameter space. 

Before we end our discussion, we would like to comment briefly on the interpretation 
of the GCE with a bino-like DM. Like the singlino-dominated DM case, a light A\ with 
mass below about 140GeV is necessary for such a work, and this A\ prefers to be singlet- 
dominated 3 . The difference is that, for the bino-like DM case, the interaction of the DM 
with A\ is relatively small and consequently the annihilation X 1 X 1 Hi,A\ can not explain 
the GCE any more due to its rather low annihilation rate. Also due to the suppressed 
interaction, rriA 1 must be closer to 2 m ^0 for the annihilation X 1 X 1 bb to account for the 
GCE, and thus the theory has to be tuned in a more elaborated way. Our sophisticated 
scan over the relevant NMSSM parameter space verified these conclusions. 

3 In the NMSSM, a light Ai with mass below about lOOGeV may have a large doublet component if the 
elements of the CP-odd Higgs mass matrix satisfy Alp ,22 3> -Mpji (see benchmark points P3 

and P4 in [122]) or M 2 P2 2 — A4 P11 ~ Afp, 12 (see the point presented in Table 2 of [123]). In either case, 
m,A 1 should be significantly smaller than m H ± to escape experimental constraints. Previous studies have 
suggested that a light doublet-dominated A\ might also explain the galactic center excess. However, due 
to the requirements on the elements this scenario occurs only in specific portions of the parameter space 
and is significantly more experimentally constrained than those we considered. In fact, in our scans for the 
GCE we did not find any parameter points with the doublet component of the light A± exceeding 0.1. In 
summary, a light doublet-dominated A\ may exist, as suggested by e.g. Ref. [123], but it is fair to say that 
without a very delicate parameter tuning, it is difficult to obtain in explaining the GCE, especially when 
one considers more constraints than previous literatures. About this conclusion, we thank the authors of 
[123] for helpful discussion. 
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A Properties of the singlet-dominated particles 


In this appendix, we present some analytic expressions for the masses and couplings of 
the singlet dominated particles, such as Xi and A\ in the NMSSM. These expressions are 
obtained by diagonalizing the mass matrices of the particles (like done in [57]), and are good 
approximations in certain cases. They are helpful in understanding the results presented 
in this work. In the following, we will follow notations and conventions consistent with 
[104] for the Z 3 NMSSM. 

A.l Neutralino masses and mixings 

In the basis ip° = (—iX\, —iX' 2 ,'ip d ,'il>u,'ips), the neutralino mass matrix is: 


/ Mi 

U 

g iv a 
V2 

0 ^ 


m 2 ^ 

g iv u 
V2 

0 


0 

“P 

-Xv u 



0 

-Xv d 

V 



2 f») 


(A.l) 


If the bino and wino fields are decoupled, the mass eigenstates of the neutralinos can be 
approximated by 


X? ~ Ni3iP° d + N U ^° U + AxsV’s, 

Xi ~ Ni + N^ s , (A.2) 


where Xi denotes the lightest neutralino with ips held as its dominant component in this 
work, and Xi represents a higgsino-like neutralino. 

In the limit of |p| » Xv, 1 » k/A and tan/3 3> 1, one can expand the neutralino 
masses and Nij by the power of Au/p = v/v s to get the following approximations: 


Ads 

Ads 

Nia 

Ads 
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cos /3 ^tan /3m ^o — p^ «-sin /3, 
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sin /3 ( p — 


m x? 


tan/3 


N 2 3 Nh \ 1/2 


1 + ) ~ Afo sa -^5gn(mjq)0(m.-.o), 


A^i 5 N 2 5 


Xv ■ R 
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V2 


1 n . ... . AT 1 An sin d . 

- 7 ^ s 9^)0(171^0), N i5 «- -= - Sgn(n)6(m^o). 

V 2 V 2 p 


(A.3) 
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In above expressions, the Sgn and 0 functions are defined by 


Sgn(x ) 


1 if x > 0, 
— 1 if x < 0, 


9( x) 


1 if x > 0, 
i if x < 0. 


(A.4) 


Likewise, one may consider the case that the wino and the singlino fields decouple. In 
this case, the mass eigenstates of the neutralinos is approximated by 


Xi ~ Nui-iXJ+Nu^ + Nu^, 

X°i « Nn(-iXi) + N i3 if° d + NuiPl (A.5) 


In the limit of tan/3 3> 1, \p\ S> g 2 V u and \/x\ S> Mi, we have the following approxima¬ 


tions: 
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A-l 
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(A.6) 


A.2 CP-odd Higgs mass matrix 

In the (A, S'/) “interaction” basis, the mass matrix for CP-odd Higgs bosons is given by 

„2 


Ml = 


m 


Mlfc S in 2/3-^) 


3/-e/i > 


AVsin2/3(^sin2/3+ |f) - 


(A.7) 


In the case of TO/i S> max(v, \A K \, \g\), k/\ <C 1 and tan/3 S> 1, the lighter CP-odd scalar 
Ai is singlet dominated with its squared mass given by 

2 9 2 ■ r) o 3 k A k I_1 

m Al ^-AK v sm 2/3--—. (A.8) 

2 A 

This approximation indicates that, without considering the radiative corrections, the singlet- 
dominated CP-odd scalar mass is determined by the parameters A, n, // as well as A K . The 
components of A\ can be written as 


Pa\.a _ 

\v 

( m A 

Pai,s / 


V 2/i 

PAi,Sj = 

f 1 + 

p2 

L4i,4 

p2 

r M C. 


si„2/3-3^) 
- 1/2 

~ 1 , 


-—sin 2/3, 

2/j 


(A.9) 


where Pa u a is the active component and Pa^.s, is the singlet component of the A\. 
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A.3 CP-even Higgs mass matrix 


In the basis (Si = cos (3h u — sin j3hd, S 2 = sin /3h u + cos /3hd, S 3 = hs ), the mass matrix el¬ 
ements for the CP-even scalars are [107] 


M\ 1 

m \ 2 

Mi 3 

M \ 2 

M 23 


M% + (rrr z — X 2 v 2 ) sin 2 2/3, 

-\( m z ~ ^v 2 ) sin4/3, 

— (M\ sin 2/3 + ) — cos 2/3, 

A (l 

m 2 z cos 2 2/3 + \ 2 v 2 sin 2 2/3, 

[ M a sin 2/3 2 k . ' 

2A ^ [‘ ~ < > *2A Sm2 ' 3 . 

2/ M Asin2/3 2 k/j 4«TC 

4 U ( /x j + T (At + ^T' 


-Akw 2 sin 2/3, 


(A.10) 

(A.ll) 

(A. 12) 
(A.13) 
(A.14) 

(A-15) 


where S 2 is nothing but the Higgs field in the SM, A4 22 is its mass at tree level without 
considering the mixing among S%, and the second term X 2 v 2 sin 2 2/3 in _A4 22 originates from 
the coupling XH U ■ HdS in the superpotential. 

The mass eigenstates Hi are defined by 


Hi = VuSi + V i 2 S 2 + V i 3 S 3 , 


(A.16) 


where V is the rotation matrix to diagonalize the mass matrix. For the S^-dominated mass 
eigenstate Hj , current Higgs data have required it to be highly SM-like, i.e. V^i, Vj )3 <C 1, 
so in the case of the hierarchy structure M.\i 3> max(Ai^ 2 , M 2 3 ), M .\ 3 -C \M \ 2 — M 2 3 1. 
If we decouple the MSSM-like heavy Higgs, Si, from the other two, the 2x2 reduced mass 
matrix in the (S 2 , S 3 ) basis is given by [107]: 

M 22 Ml 2 3 

M 33 — A ie^ A s ™ 2 4/3 — K ^ 2 ° cos 2 2/3 — Xkv 2 cos 2 2/3 sin 2/3 

The (2,2) element of the reduced (2 x 2) matrix, which in the limit of zero-mixing with the 
other Higgs should give singlet scalar mass in the Z 3 NMSSM, is given by: 




( 2 ? 2 ) — -y- ( + 


A 


4k/A A 2 v 2 m\ 


+ 


4 /x 2 


(l — cos 2 2/3) sin 2 2/3 


■ K cos 2 2/3 — -kXv 2 (2 cos 2 2/3 + l) sin2/3. 
m , 2 


(A.18) 


Setting M 23 ~ 0, i.e. m 2 A = si ^ 2 2 /3 ~ JX sin 2/3); and taking tan/3 1, we have 

Ml s ,( 2,2).f(^ + ^) (A. 19) 


This approximation indicates again that, without considering the radiative corrections, the 
singlet-dominated CP-even scalar mass is determined by the parameters A, k, /x and A K . 
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A.4 Some properties of the singlet-dominated particles 


With the assumptions that M A S> max(\g\,\A K \), \/jl\ 3> Xv, tan/3 1 and k/X <C 1, 
one can approximate the masses and couplings of the singlet dominated particles, such 
as Xi and A\ , by simple analytic expressions [57]. In the following, we list some of the 
coupling expressions used in our discussion, which are denoted by Cxyz hereafter. These 
expressions are actually expand the corresponding exact ones by the power of Xv/g. 

im\, tan f5 inn, Xv 

Ca ^ = ^T Pma "~VI vj* 

_ j -iV2n(l + 2^) for j = 1, 

AiXiXj ~ | “t7T ^gn(m^q)9(mg>) for higgsino — like Xj, 

Cjj.£ 0^,0 = —iC Al ^o^o if Hj is singlet dominated. (A.20) 

Likewise, if Hi is the SM-like Higgs boson, we have 

C HiX.iXj 


T,Sgn{m^)6(m^o) 


for j = 1, 

for higgsino — like Xj, 


(A.21) 
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